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Abstract 

In  this  paper  we  address  theoretical  issues  related  to  a  class  of  nonlin¬ 
ear  parabolic  systems  used  in  toxicokinetic  modeling.  We  establish  the  well- 
posedness  of  these  systems  in  a  weak  or  variational  setting  by  proving  the 
existence  of  a  unique  solution  which  depends  continuously  on  the  data.  These 
results  are  then  applied  specifically  to  the  model  equations  for  a  physiologically 
based  pharmacokinetic  model  for  the  systemic  transport  of  trichloroethylene 
(TCE),  which  is  the  motivating  example  for  our  efforts  here. 

In  addition,  we  present  theoretical  results  for  the  parameter  estimation 
problem  associated  with  the  general  class  of  abstract  systems  and  the  example 
TCE  model.  These  results  further  establish  the  convergence  of  the  Galerkin 
approximations  used  in  the  numerical  scheme  for  the  TCE  model. 


1  Introduction 


In  this  paper  we  present  theoretical  results  related  to  a  class  of  physiologically  based 

pharmacokinetic  (PBPK)  models.  PBPK  models  are  widely  used  in  the  field  of 

toxicology  to  describe  the  uptake,  systemic  transport  and  elimination  of  compounds 

within  humans  and  animals.  These  compartmental  models  lead  to  a  system  of 
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coupled  algebraic  and  differential  equations,  where  each  equation  represents  the 
transport  of  the  compound  through  a  given  organ  or  tissue.  The  standard  PBPK 
models,  known  as  the  perfusion-limited  and  diffusion-limited  models,  are  ODE- 
based  models  that  are  derived  using  simple  mass  balance  principles  and  assumptions 
of  uniformity  and  rapid  equilibrium  [20]. 

In  some  cases,  the  assumption  of  a  “well-mixed”  compartment  is  not  appro¬ 
priate  for  tissues  with  heterogeneous  physiological  properties.  For  example,  large 
degrees  of  spatial  variations  in  cell  size,  lipid  distribution,  and  membrane  perme¬ 
ability  within  adipose  (fat)  tissue  may  lead  to  spatially-dependent  concentrations 
of  highly  lipophilic  (fat-soluble)  compounds  there.  To  more  accurately  capture  the 
effects  of  these  tissue  heterogeneities,  specialized  PDE-based  compartmental  mod¬ 
els  can  be  developed  by  incorporating  the  specific  physiology  of  the  tissue  into  the 
model  assumptions.  The  resulting  compartment  is  then  coupled  with  the  standard 
PBPK  compartments  for  the  remaining  organs  and  tissues,  yielding  a  hybrid  PBPK 
whole-body  model. 

Examples  of  hybrid  PBPK  models  include  a  CFD-PBPK  model  for  nasal  gas 
and  vapor  uptake  [12],  a  distributed  parameter-PBPK  model  for  inhalation  and 
dermal  absorption  [26],  and  a  hybrid  PBPK  model  for  the  systemic  transport  of 
trichloroethylene  (TCE).  The  TCE  model,  which  motivates  our  efforts  here,  includes 
a  specialized  axial  dispersion-based  model  for  the  transport  of  TCE  in  the  adipose 
tissue,  where  TCE  is  known  to  accumulate  [29]. 

TCE  is  a  highly  lipophilic  solvent  that  has  been  used  widely  in  industry  as  a 
metal  degreasing  agent.  Now  a  common  soil  and  groundwater  contaminant,  TCE 
can  be  found  at  Superfund  sites  and  Department  of  Defense  facilities  across  the 
United  States  [28].  Humans  come  into  contact  with  TCE  most  often  by  inhaling 
TCE  vapor  in  an  industrial  setting,  or  more  commonly  by  drinking  contaminated 
water. 

TCE  and  several  of  its  metabolites  are  known  to  produce  toxic  effects  in  humans. 
Acute  exposure  to  TCE  impacts  the  central  nervous  system,  causing  symptoms  such 
as  fatigue,  headaches,  dizziness  and  drowsiness  [16].  Several  of  its  metabolites  are 
suspected  of  causing  toxicity  in  animals  and/or  humans  [10,  21],  including  liver, 
lung  and  kidney  tumors  [9,  10,  11,  17],  as  well  as  developmental  defects  [8,  27]. 

A  key  property  of  TCE  that  greatly  affects  its  pharmacokinetics  is  its  high  degree 
of  lipophilicity,  which  leads  to  the  accumulation  of  TCE  in  the  adipose  tissue.  Nu¬ 
merous  studies  have  demonstrated  that  adipose  tissue  is  highly  heterogeneous,  with 
wide  variations  in  fat  cell  size,  metabolic  activity,  blood  flow  rates,  and  membrane 
permeability  [13,  14,  18,  25].  Moreover,  the  amount  of  lipid  in  each  adipocyte  is 
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known  to  vary,  creating  an  uneven  distribution  of  lipids  across  the  tissue  [18].  Since 
TCE  and  other  compounds  of  interest  are  highly  lipophilic,  this  uneven  distribu¬ 
tion  of  lipids  may  lead  to  an  uneven  distribution  of  the  compounds  themselves. 
This  would  then  result  in  a  spatially  dependent  concentration  function  for  such  a 
compound  inside  the  adipose  tissue. 

In  Section  2  we  present  an  overview  of  the  TCE  PBPK-hybrid  model,  including 
the  state  space  settings  and  weak  formulation.  This  model  is  the  motivation  for  the 
theoretical  results  we  develop  here  for  a  general  abstract  class  of  nonlinear  parabolic 
systems  that  may  arise  in  toxicokinetic  modeling.  In  Section  3  we  establish  the 
existence  of  unique  weak  solutions  for  the  general  class  of  problems  and  the  TCE 
PBPK-hybrid  model.  Finally,  in  Section  4  we  address  theoretical  issues  related  to 
the  associated  parameter  estimation  problem.  The  results  presented  in  this  paper 
also  provide  a  framework  for  a  numerical  scheme  that  we  implement  computationally 
in  [22]  for  model  simulation  and  parameter  identification.  Moreover,  we  establish 
the  theoretical  convergence  of  numerical  solutions  to  the  solution  for  the  infinite¬ 
dimensional  system  of  model  equations. 


2  The  TCE  PBPK-hybrid  model 

In  this  section  we  provide  an  overview  of  the  TCE  PBPK-hybrid  model,  which  is 
developed  and  discussed  in  full  detail  in  [2,  22],  The  compartmental  model  for  the 
adipose  tissue  is  based  on  the  axial  dispersion  model  of  Roberts  and  Rowland  [24] 
for  the  transport  of  solutes  in  the  liver,  and  is  adapted  for  the  specific  physiology 
of  adipose  tissue.  This  aggregate  model  uses  a  representative  “cell”  to  capture 
the  transport  behavior  in  a  collection  of  many  similar  “cells”  that  have  varying 
properties. 

In  this  particular  case,  the  representative  unit  includes  a  single  adipocyte  (fat 
cell)  with  an  adjoining  capillary,  both  immersed  in  the  surrounding  interstitial  fluid. 
As  discussed  in  [2,  22] ,  the  adipocyte  is  represented  by  a  sphere  centered  at  the  origin 
with  radius  r\ ,  and  the  capillary  is  a  cylindrical  tube  with  circular  cross-section  and 
central  axis  r  =  r-2+a,  8  =  8q,  E\  <  cj>  <  e2  in  spherical  coordinates.  The  interstitial 
fluid  fills  in  the  space  surrounding  the  other  two  regions. 

We  assume  that  blood  and  TCE  enter  the  representative  adipose  capillary  region 
from  the  arterial  blood  system.  The  transport  within  the  capillary  is  described  by 
a  one-dimensional  convection-dispersion  equation,  while  the  exchange  between  the 
capillary  and  the  other  two  adipose  regions  (adipocyte  and  interstitial  space)  is 
modeled  using  Fick’s  first  law  of  diffusion.  Moreover,  we  use  flux  balance  to  derive 
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boundary  conditions  that  connect  the  capillary  region  to  the  arterial  and  venous 
blood  systems.  Within  each  of  the  adipocyte  and  interstitial  regions,  we  have  a 
two-dimensional  diffusion  equation,  as  well  as  the  terms  representing  the  passive 
transport  between  all  three  regions.  The  boundary  conditions  in  the  adipocyte  and 
interstitial  space  are  standard  periodic  and  finiteness  conditions  that  are  appropriate 
for  diffusion  on  a  spherical  domain. 

The  other  tissue  compartments  used  in  the  PBPK-hybrid  model  correspond 
to  the  brain,  kidney,  liver,  muscle  and  remaining  tissues,  and  each  is  modeled 
as  a  perfusion-limited  compartment.  We  assume  uptake  via  inhalation  in  the 
lungs,  which  is  modeled  using  a  steady-state  assumption.  Moreover,  we  include 
a  Michaelis-Menten  term  in  the  liver  to  represent  the  metabolism  of  TCE  there. 
See  [2,  22]  for  a  detailed  derivation  and  discussion  of  the  model  equations. 

2.1  State  space  setting  and  weak  formulation 

Here  we  present  the  weak  or  variational  formulation  for  the  TCE  PBPK-hybrid 
model.  First  we  consider  the  three  regions  of  the  adipose  tissue  compartment.  The 
domain  for  the  capillary  region  of  the  adipose  tissue  is  given  by  the  arc 

=  {{r,9,(f>)  :  r  =  n,9  =  0o,£i  <  <f>  <  ir  -  s2} 

in  spherical  coordinates,  and  the  domain  for  the  adipocyte  and  interstitial  regions 
is  the  spherical  shell 

Qia  =  {(r,  9,(j>)  :  r  =  r\ ,  0  <  9  <  2ir,  0  <  <p  <  7 r}. 

Note  that  both  r  and  9  are  constant  in  the  domain  fls,  so  that  we  may  write 
u  =  u((j>)  for  functions  u  G  Qb-  Similarly,  since  r  is  constant  in  Qja  we  write 
u  =  u(9,  <j>)  for  u  G  &ia-  In  both  cases,  r  =  r\  is  understood. 

We  choose  the  following  state  spaces  for  the  regions  of  the  adipose  dispersion 
model:  Hg  =  H1(fls)  in  the  capillary  region  and 

H}a  =  Hper(flIA)  =  {u  G  H^CIia)  :  u(9,  <j>)  =  u(9  +  2tt,^)} 

in  the  interstitial  and  adipocyte  regions.  The  norms  and  inner  products  on  these 
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spaces  are  defined  by 
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Note  that  in  spherical  coordinates,  for  u  G  H1(flB)  we  have 


U  =  \Jh'  r  sin  0  ~d6  '  r  "d<f>  J  =  V  ’ 

and  for  u  G  4)  we  have 
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We  now  define  the  state  space  V  =  x  H}a  x  H}a  x  K6  and  the  space 
Ti  =  L2(ClB)  x  L2(CliA)  x  L2(£liA)  x  M6,  with  corresponding  norms  and  inner 
product 
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Let  T  denote  the  boundary  of  the  domain  .  Since  0/^  is  defined  as  a  spherical 
shell  in  spherical  coordinates,  we  can  choose  the  boundary  T  to  be  any  fixed  value 
of  6.  Therefore  we  define  T  =  {{r,  9,  d)  ■  r  =  r±,  9  =  9o,  0  <  <  7r}.  Moreover,  we 

define  the  trace  operator  Tgo  :  H1( fljA)  — >  L2(T)  so  that 

Tlu=u\eS  0  (1) 

for  any  u  G  4),  where 

\u\l2<d  =  f  \u(<f>)\2  sin  <f>  # 

JO 
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for  all  u  G  T2(T),  and  where  \u\hi^qia^  =  \u\Hi  . 

It  follows  from  the  Trace  Theorem  (see  Theorem  6  in  [23],  p.  240)  that  there 
exists  K\  >  0  such  that  for  every  u  G  H1(Qja)  we  have 
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Here  we  are  seeking  a  solution  y(t)  G  V  which  satisfies  the  initial  condition 
J/(0)  =  yo  G  7i  and 
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The  vector  y  represents  TCE  concentrations  (in  mg/liter)  in  each  of  the  tissue 
compartments  and  subcompartments.  Specifically,  we  have 

y(t)  =  [ CB(t ),  Ci(t),  CA(t),  cv(t ),  Cbr(t),  Ck(t ),  Cm(t),  Ct(t),  Cl(t)]T 

where  CbH),  Ci(t )  and  CA(t)  represent  TCE  concentrations  in  the  capillary,  inter¬ 
stitial  and  adipocyte  regions  of  the  adipose  tissue  respectively,  and  Cv{t),  Cbr{t), 
Ck(t),  C  ,„(!■),  Ci(t )  and  Ct{t )  denote  TCE  concentrations  in  the  venous  blood, 
brain,  kidney,  muscle,  liver  and  remaining  tissue  compartments,  respectively.  The 
function  Cc(t)  that  appears  in  the  forcing  function  (5)  is  the  concentration  of  TCE 
in  the  chamber  air,  and  is  a  specified  quantity.  See  [2,  22]  for  details. 

Here  we  use  the  notation  y  =  ||,  and  {■,  *)y,v  denotes  the  usual  duality  product 
described  in  the  next  section.  The  system  (2)  belongs  to  a  general  class  of  abstract 
problems  that  we  consider  in  the  next  section,  where  we  address  issues  of  well- 
posedness. 


3  Well-posedness  of  solutions 

First  we  present  well-posedness  results  for  a  general  class  of  abstract  nonlinear 
parabolic  equations.  The  theoretical  issues  that  we  address  here  establish  the  ex¬ 
istence,  uniqueness  and  continuous  dependence  of  solutions  for  the  TCE  PBPK- 
hybrid  model,  as  well  as  for  a  class  of  more  general  equations  that  may  arise  in 
more  complicated  models  of  xenobiotic  transport  inside  the  body.  Moreover,  these 
theoretical  efforts  lay  the  groundwork  for  computational  methods  used  in  model 
verification,  simulation  and  parameter  estimation. 

The  main  result  that  we  establish  here  regarding  well-posedness  of  solutions  is 
based  on  ideas  presented  in  [5]  and  [1] .  Banks  and  Musante  [5]  proved  well-posedness 
results  for  two  classes  of  abstract  nonlinear  parabolic  systems.  The  first  class  they 
addressed  requires  the  nonlinearity  to  satisfy  a  particular  convexity  condition.  The 
second  class  replaces  the  convexity  condition  with  an  assumption  that  the  forcing 
function  has  additional  regularity  with  respect  to  time  or  space. 

The  TCE  PBPK-hybrid  model  has  a  simple  Michaelis-Menten  nonlinearity  in 
the  liver  equation,  and  it  can  be  shown  that  this  system  of  model  equations  fits 
into  the  first  of  these  abstract  classes  of  problems.  Although  this  happens  to  be 
the  case  for  TCE,  which  has  a  simple  metabolic  mechanism,  there  may  be  other 
lipophilic  compounds  which  have  more  complicated  kinetics  (e.g.,  protein  binding) 
that  would  not  satisfy  the  convexity  condition. 

The  forcing  function  in  our  TCE  model  is  related  to  the  inhalation  of  TCE  va- 
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por  in  the  air,  and  does  not  have  the  necessary  regularity  required  for  the  model 
to  fit  into  the  second  class  of  abstract  problems  discussed  by  Banks  and  Musante. 
Transport  models  for  other  inhaled  compounds  would  have  a  similar  forcing  func¬ 
tion  that  also  may  lack  the  required  regularity.  Therefore,  a  model  for  an  inhaled 
lipophilic  compound  with  complicated  nonlinear  kinetics  would  likely  fall  outside 
both  of  these  abstract  classes  of  problems. 

Using  ideas  in  [1],  we  can  improve  upon  the  results  of  Banks  and  Musante  by 
achieving  well-posedness  for  a  more  general  class  of  abstract  nonlinear  parabolic 
equations.  Ackleli,  Banks  and  Pinter  [1]  proved  the  well-posedness  of  second  order 
elastomer  problems  in  which  the  nonlinearities  satisfy  local  Lipschitz  and  affine 
domination  properties.  Our  result  is  based  on  their  work,  and  will  establish  well- 
posedness  for  a  wider  class  of  inhaled  xenobiotic  transport  models. 


3.1  Well-posedness  results  for  a  class  of  abstract  nonlinear 
parabolic  systems 


Consider  the  system 

y(t)  +  Ay(t)  +  g(y(t))  =  f(t)  in  V*  (7) 

2/(0)  =  2/o  (8) 

for  t  E  (0,T)  with  T  <  oo.  We  assume  that  V,  %  and  V*  are  separable  real  Hilbert 
spaces  that  form  a  Gelfand  triple  [30]  which  satisfies 

V  ^  n  ~  n*  V*. 

We  denote  the  inner  product  in  %  by  {•,•),  and  the  norms  in  each  of  the  spaces  V, 
1~L  and  V*  are  denoted  by  |-|v,  |-|  and  |-|v.  respectively.  Moreover,  we  assume  that 
the  embedding  V  T-L  is  dense  and  continuous,  with 

\ip\  <  k  \ip\v  for  all  ip  €  V.  (9) 

The  duality  product  (•,  -)v*,v  is  the  extension  by  continuity  of  the  inner  product  in 
%  from  Ti  x  V  to  V*  x  V  (see  [30]  for  a  complete  discussion). 

We  assume  the  operator  A  is  defined  in  terms  of  a  given  sesquilinear  form  a  : 
V  x  V  — >  R.  That  is,  we  define  A  :  V  — >  V*  by  (Au,v)v,v  =  <r(u, i>),  and  it  follows 
under  the  assumptions  on  a  below  that  A  G  £(V,  V*). 

We  make  the  following  standing  assumptions  on  a,  /  and  g: 


(Al)  The  sesquilinear  form  a  is  bounded  in  V,  i.e.,  there  exists  Ci  >  0  such  that 


\a(u,v)\  <  Gi  |w|v  Hv 
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(10) 


for  all  ii,  v  G  V. 


(A2)  The  sesquilinear  form  a  is  elliptic  on  V.  That  is,  there  exist  k\  >  0  and  A  >  0 
such  that 

o(u,u )  >  k\  \u\~v  —  A  \u\~  for  all  u  G  V.  (11) 

(A3)  The  forcing  function  /  satisfies 

/ei2((0,T);V*).  (12) 

(A4)  The  nonlinear  function  g  :  72  — )■  72  satisfies  the  following  local  Lipschitz 
condition:  let  Br( 0)  =  {u  G  72  :  |w|  <  r}  denote  the  ball  of  radius  r  centered 
around  the  origin  in  Li.  Then  given  r  >  0,  there  exists  Lsr  >  0  such  that 


| g(u)  —  g(v) |  <  Lsr  |u  —  u|  for  all  u,v  G  Br( 0).  (13) 

(A5)  There  exist  positive  constants  Cn  and  C3  such  that 

\g(u)\  <  Co  |m|  +  C3  for  all  u  G  H.  (14) 

We  say  that  y  G  T2((0,T);  V)  is  a  solution  of  (7),  (8)  if  it  satisfies 

(yW^)v,v  +p(v(t),$)  +  (g(y(t)),^)  =  (/(*},  ^)v*,v  (15) 

2/(0)  =  i/o  (16) 


for  all  c  G  V. 

Theorem  3.1  (Local  existence)  Under  Assumptions  (Al)  -  ( A 4)  and  for  any  yo  G 
72,  the  system  (15)  -  (16)  has  a  unique  solution  y  G  LodO.T);  V)  fl  C.([0,T];  72)  on 
some  interval  [0,2*]. 

Proof:  Let  V  denote  the  Hilbert  space  radial  retraction  onto  the  ball  in  72  with 
radius  1  centered  at  y 0,  and  define  the  nonlinear  function  g  :  72  — >  72  by 

g(u)  =  g(Vu), 


Then 

I  g(u)  -  g(v)  |  <  LB(  |„0|+1)  1  Vu  -  Vv\  <L\u-  v\  (17) 

for  all  u,v  G  72,  where  L  =  2  (see  Lemma  2.1,  p.  230  in  [4]).  Moreover, 

\g(u)\  <  L\u\  + c  (18) 

where  c  =  |<?(0)|  >  0  since 

\g(u)  \  -  | <)(0)|  <  | g(u)  -  3(0)|  <  L  |u| . 
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Now  consider  the  problem 

y{t)  +  Ay(t)  +  g(y(t ))  =  f(t)  in  V*  (19) 

2/(0)  =  yo,  (20) 

which  is  equivalent  to 

+  p(v(t),4)  +  {g(y(t)),^)  =  {f(t)i.$>)v,v  (21) 

2/(0)  =  yo  (22) 


for  all  ip  e  V. 

Suppose  that  there  exists  a  solution  y(t)  to  (21)  -  (22)  on  an  interval  [0 ,T], 
Here  we  establish  a  priori  bounds  for  this  solution  by  setting  ip  =  y(t )  in  (21)  and 
integrating  from  0  to  t  for  some  t  £  [0,T]: 

\  1 2/(1) |2  +  [&{y(s),y(s))  +  (g(y(s)),y(s))]ds  =  ^  |y(0)|2  +  (f(s),y(s))v,,vds. 

Moreover,  using  (11),  (18),  (22)  and  standard  inequalities  we  obtain 

|'2/(1)|2  +ki  f  [y(s)\2vds  <  |'2/o|2  +  t-  /  \f(s)\2v,  ds  +  c2T  +  LCT 

Jo  K\  Jo 

=  C(ki,  f,  yo,c,T,  L,C).  (23) 


This  establishes  an  a  priori  bound  for  any  solution  of  the  system  (21)  -  (22). 

Next  we  define  the  “Galerkin”  approximations  for  (19)  by 

N 

yN  (i)  = 

k= i 

where  {/•.<},)  ,  C  V  is  a  linearly  independent  total  subset  of  V,  and  where 
{cf  (!)}0=1  are  chosen  so  that  yN (t)  is  the  unique  solution  of 

(yN(tUj)+a(yN(t),rPj)  +  ( g(yN (t ))&)  =  </(/).  r,)y  .v  (24) 

for  j  =  1, . . . ,  N  with  c^P (0)  =  q^..  Since  {ipk}kLi  is  a  total  subset  of  V,  we  can  and 
do  choose  the  constants  {q^,}  so  that 

N  N 

Vo  =  Jim  H  yo  =  Jim  c0k^k  in  V.  (25) 

N-^oo  z — '  N— >oo  z — ' 

A:— 1  k=l 

Now  we  multiply  (24)  by  Cj1  (t)  and  sum  over  j  =  1, . . . ,  N: 

\jt  +  o(yN(t),  yN(t))  +  (g(yN(t)UN(t))  =  (M,yN(t))v.,v. 


10 


Note  that  \'!Jo  \v  is  uniformly  bounded  by  some  constant  A'o  since  converges 
strongly  to  yo  as  in  (25),  and  therefore  we  can  use  arguments  similar  to  those  above 
to  obtain  the  uniform  bound 

|  yN(t)\2+h  f  \yN{s)^ds  <c{k,J,K^T,L,C)  (26) 

Jo 

for  some  c  which  is  independent  of  N. 

Next  we  turn  to  the  convergence  of  the  Galerkin  approximations.  Note  that  (26) 
implies  that  {i/A  }  is  bounded  in  L'2((0,T);V).  Therefore  there  exists 
y  G  L'2((0,T);  V)  and  a  subsequence  of  {yA  }  (which  we  will  denote  by  {«/A  })  such 
that 

//'  m  T2((0,T);V).  (27) 

It  follows  that  yN  y  in  L2((0,T);  H),  which  further  implies 

yN  (t)  y(t)  in  T-L  for  almost  every  t  €  [0,T].  (28) 

Moreover,  (18)  and  (26)  imply  that  {g(yN )}  is  bounded  in  I/2((0,T);  A),  so  that 
there  exists 

h  G  L2((0,T);?{)  such  that 

g(yN)^h  in  L2((0,T);A).  (29) 

Next  we  establish  that  y  is  a  solution  to  (21),  (22).  Let  Pm  be  the  class  of 
functions  G  L2((0,T);  V)  that  can  be  represented  in  the  form 

M 

v(t)  =  ^2ak{t)ipk, 

k=l 

oo 

where  a&  G  C1[0,T].  Moreover,  we  define  P  =  Pm ,  which  implies  that  P  is 

M= 1 

dense  in  L2((0,  T);  V). 

We  multiply  (24)  by  aj(t),  then  sum  from  1  to  M  and  integrate  over  [0 ,t], 
integrating  by  parts  in  the  first  term.  Then  for  fixed  G  Pm  for  M  <  N,  we  take 
the  limit  as  N  — >  oo  to  obtain 

/  [-(y(s),v(s))  +  a(y(s),ii(s ))  +  (h(s),y(s))]ds  +  (y(t),i](t)) 

Jo 

=  <2/o,77(0))  +  /*  {f(s),r](s))vyds  (30) 

Jo 

for  almost  every  t  G  [0,  T],  where  we  have  used  the  convergences  (27),  (28)  and  (29). 
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Note  that  the  M  we  chose  above  was  arbitrary,  and  since  {ipj}  is  total  in  V  it 
follows  that 

,T  ,T 


/  d{s){-y{s},ip)ds  +  / 


a(s)  [ a(y(s),ip )  +  {h(s),ip)  -  (f{s), .ip)v,v]  ds  =  0 


Jo  Jo 

for  every  ip  G  V.  This  equation  holds  for  every  a  G  Cq°[0,T],  which  further  implies 
that  y  G  L2((0,T);  V)*  ~  La((0,T);  V*)  and 


(: y(t),ip)v\v  +  +  < h{t),ip )  =  {f(t),tp)v,v  (31) 

for  each  ip  G  V  and  for  almost  every  t  G  [0,  Tj.  Moreover,  since  y  G  T2((0,  T);  V)  and 
y  G  L2((0,T);  V*),  it  follows  from  Theorem  3.1  in  [19,  p.  19]  that  y  G  C([0,  T];  %). 

To  show  that  //(())  =  yo,  consider  (30)  with  y  =  yj  =  a(t)ipj,  a  G  C'1[0,T]  and 
a(0)  0,  integrating  by  parts  in  the  first  term: 

/  [(y(s),/?(s))v*,v  +  <r(y(s),il(s))  +  (/i(s),'/?(s))]ds  +  (y(t),y(t)) 

JO 

+  <-y(s),r?(s))|^o  =  <2/o,77(0))  +  /  {f(s),y(s))v%vds.{ 32) 

Jo 

It  follows  from  (31)  and  (32)  that  y(0)  =  yo- 

Finally,  we  show  that  h  =  #(■(/).  If  we  let  2A  (t)  =  yA  (i)  —  j/(i)  G  L2((0,T);  V)  fl 
C([0,T];  H),  and  combine  (24)  with  ipj  =  yN  ( t )  and  (31)  with  ip  =  j/(i),  integrating 
from  0  to  I  and  adding,  we  obtain 

|zJV(*)T+2  [  (r(zN(s)i.zN(s))ds  =  |i/JV(0)  -  y0\~  +  2  /  {f(s),zN(s))v,,vds 

Jo  Jo 

+  4/  {f(s),y(s))v*,vds  -  2  f  (g(yN(s))  -  h(s),  zN(s))ds 

Jo  Jo 

-  2  f \h(s),yN(s))ds  -  2  [\g{yN{s)),y{s))ds  -  2 (y(t),yN(t)) 

Jo  Jo 

+  2(2/0, 2/^(0))  -  2  f  a(yN(s),y(s))ds  -  2  f  a(y(s),yN (s))ds.  (33) 

«/o  </o 


Note  that 

2  f \g(yN(s))-h(s),zN(s))ds 
Jo 


=  2  [\g(yN(s))-g(y(s)),zN(s))ds 

Jo 

+  2/  (g{y(s))  —  h(s) ,  zN  (s))  ds , 


J  0 

and  by  applying  (17)  and  (11)  to  (33)  we  obtain 

|zjv(I)|2+2  fh  \zN(s)\lds  <  \yN(0)-yo\2 +  2(L  +  X)  f 

Jo  Jo 

+  |xJV(f)|  +  |iaJV(i)|  , 


^^(s)  |  ds 

(34) 
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where 


XN{t)  =  2  —{y{t),yN{t))—  f  a(;y(s),yN (s))ds  —  f  a(yN  (s),y(s))ds 

Jo  Jo 

+  {yo,yN{0))~  f  (g{yN{s)),y{s))d,s  -  f  (h{s),yN{s))ds 

Jo  Jo 

+  2  I  (f(s),y(s))v,vds  , 

Jo 

yN(t)  =  2  I  (g(y{s))  -  h(s),zN(s))ds  +  2  f  (f(s),zN(s))v*,vds. 

Jo  Jo 


Note  that  (25)  and  (27)  imply  that  |;i/A(0)  —  yo\~  — »  0  and  |l'JV(i)|  — »  0  respec¬ 
tively.  It  also  can  be  shown  that  |XA(i)|  — »  0  by  taking  the  limit  as  N  — y  oo  and 
using  the  convergences  (28),  (27),  (25)  and  (29): 


XN{t)  ->•  2  -  [y{t)\2  -  2  f  a(y(s),y(s))ds  +  \y0\2 

Jo 

-2  /  (h(s),y(s))ds  +2  [  ( f(s),y(s))vtVds 


Jo  Jo  J 

which  is  two  times  the  integrated  form  of  (31)  with  ib  =  y.  Therefore  |A"Ar(^|  — y  0 
and  (34)  becomes 


\z*  (t) \~  +  2  /  k\  | z^  (s)|“  ds  <  L  f  | ^^(s) | -  ds  +  WN ( t ), 

Jo  Jo 

where  L  =  2 (L  +  A),  and  where  WN ( t )  =  \yN (0)  —  yo\"  +  |dVJV(i) |  +  |yJV(i)|  — >  0 
as  N  — >  oo.  Using  GronwalPs  inequality  we  obtain 

\zN(t)\2  <  WN(t)  +  I  LWN{s)eLt-Lsds  ->  0, 

Jo 

which  implies  that  \z^  (i)|“  — >  0  for  t  G  [0,T]  and  hence  yN(t)  — »  y(t )  in  U  for  each 
t.  It  follows  that 

/  | g(y*  (s))  -  g(y(s))\~  ds  <  L2  f  \yN(s)  -  y{s)\~  ds  ->■  0 
Jo  Jo 


so  that  g(yN )  — »  g(y)  in  Z/2((0,T);%).  Since  (29)  implies  that  g(yN)  —^h  in 
we  have  h  =  g(y)  in  Z,2((0,  T);  H)  and  therefore 


(h,ip)  =  (g(y),il>) 

for  every  ip  G  V.  This  establishes  that  y  is  a  solution  of  (21),  (22).  The  uniqueness 
of  solutions  can  be  shown  using  standard  arguments  (see  [3,  5]),  and  is  presented 
in  detail  in  [22]. 
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Finally,  we  show  the  existence  of  a  solution  for  the  original  problem  (15),  (16) 
on  some  interval  in  time.  Since  the  solution  y  of  (21),  (22)  satisfies  y  G  C*([0,  Tf|  72), 
there  exists  t*  £  [0,T]  such  that 


I y(t)  —  2/o  |  <  1 

for  all  t  G  [0,t*].  Therefore  y(t)  G  Bi(y0)  for  all  t  G  [0,2*],  which  further  implies 
that  g(y(t ))  =  g{Vy(t))  =  g{y(t))  for  all  t  G  [0,2*].  It  follows  that  y  is  a  solution 
of  (15),  (16)  on  the  interval  [0,2*].  The  uniqueness  of  the  solution  y  can  be  shown 
using  standard  arguments.  This  completes  the  proof  of  the  theorem. 

The  existence  of  a  unique  global  solution  to  (15),  (16)  now  can  be  shown  using 
the  additional  assumption  (A5). 

Theorem  3.2  (Global  existence )  Under  Assumptions  (Al)  -  (A5)  with  yo  G  72, 
the  system  (15),  (16)  has  a  unique  global  solution  y  G  Lo((0,T);V)  H  C([0,T];72) 
which  depends  continuously  on  the  data  {yo,f)- 

Proof:  Let  yN(t)  =  J2k=i  ck  (Oi’k,  where  {cj^ ( t )}  are  chosen  so  that  yN (2)  is  a 
solution  of  the  finite  dimensional  system  of  equations  derived  from  (15),  (16)  with 
ipk  as  test  functions.  Moreover,  we  choose  the  constants  {c^}  so  that 

N 

yo  =  lim  yo  =  lim  r'k '  '*•  in  V.  (35) 

N— >oo  z — '  A'— s-oo 

k=l 

As  before  we  obtain 

[yN{t)\~ +2  f  a(yN{s),yN{s))ds  +  2  /  {g(yN (s)),yN (s))ds 

Jo  Jo 

=  2/  {f(s),yN(s))v%vds +\yN(0)\\ 

Jo 

and  using  arguments  similar  to  those  for  Theorem  3.1  we  arrive  at  the  uniform 
bounds 

\yN(V\2  <C(yo,ki,C3,T,C2,X,f,C0)  (36) 

and 

| yN {t)\'+h  [  | yA"  (s)  |  ^  ds  <  Cx 
Jo 

for  t  G  [0,T],  where  T  is  arbitrary. 

This  implies  that  there  exists  y  G  L-2((0,T);V )  and  a  subsequence  of  {;yAr}  such 
that 

VN^  y  in  Lo((0,T);V) 
yA  ( t )  y(t)  a.e.  in  7 i 
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g(yN)^h  in  L2((0,T); 70- 

In  addition,  we  use  the  weak  lower  semicontinuity  of  norms  to  obtain 

\y(t)\ 2  <  liminf  \yN  (*)|“  <  C  (37) 

iv— >-00 

for  t  G  [0 ,T],  Therefore  y(f)  G  5^r(0)  in  H  for  0  <  t  <  T,  and  the  arguments  of 
Theorem  3.1  imply  that  y  is  the  unique  global  solution  to  (15),  (16). 

To  prove  continuous  dependence,  we  consider  solutions  y  and  y  in  L2((0,T);  V) 
to  (15),  (16)  with  initial  conditions  yo  and  y0  and  forcing  functions  /  and  /,  respec¬ 
tively.  Note  that  since  y  G  L2((0,T);V)  satisfies  (15)  with  y( 0)  =  yo,  we  can  use 
arguments  similar  to  those  used  earlier  in  this  proof  for  yN  (t)  to  obtain  the  bounds 

\y(t)\ 2  <  C(y0,k1,C3,T,C,,f,\)=C, 

|y(t)|2  <  C(y0,k1,C3,T,C,,f,X)=C3. 


Therefore  y(t)  G  anc^  e  ^ which  imply  that 

\g(y(t))  -g(y(t))  \  <  lBiS4  \ y(t)  -y(t)\, 

where  C'4  =  max{  \ZCo,  x/c^}- 


(38) 


Now  consider  w(t)  =  y(t)  —  y(t),  which  satisfies  w(0)  =  yo  —  jjo  and 

(w{t),ip)v*y  +  (39) 

for  all  ip  G  V,  where  A g(t)  =  g{y{t))  —  g{y{t))  and  /(f)  =  /(f)  —/(f)-  Note  that  (39) 
is  also  satisfied  for  all  ip  G  L2((0,T);  V). 


For  fixed  f  G  [0,T]  we  choose  ip  =  w(t )  G  V  in  (39)  and  integrate  from  0  to  f: 


[(w(s),w(s))v*,v  +<r (w(s),w(s))  +  (A g(s),w(s))]ds  =  /  (f(s),w(s))v*yVds. 


We  then  apply  the  conditions  (11),  (38)  and  GronwalPs  inequality  to  obtain 


|w(f)|  <  e  c4 


H0)\ '2  +  y  1 

h  Jo 


f(s) 


ds 


for  all  f  G  [0,T],  Therefore  we  have  shown  continuous  dependence  in  G([0,T];7f) 
on  the  initial  data  (yo,  /).  This  completes  the  proof  of  Theorem  3.2. 


3.2  Well-posedness  for  the  TCE  model 

In  this  section  we  verify  that  the  TCE  model  fits  into  the  class  of  abstract  nonlinear 
parabolic  systems  discussed  in  the  previous  section.  This  will  establish  the  well- 
posedness  of  the  TCE  model. 
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Consider  the  weak  formulation  (2)  of  the  TCE  model  with  the  definitions  of  V, 
H  and  the  corresponding  norms  and  inner  products  as  outlined  in  Section  2.1.  The 
spaces  V  and  7i  form  a  Gelfand  triple  as  discussed  in  Section  3.1  with  the  duality 
product  (v)y,v-  It  can  be  shown  using  standard  arguments  that  Assumptions 
(Al)  -  (A5)  are  satisfied  for  the  TCE  PBPK-hybrid  weak  formulation  (see  [22]  for 
detailed  proofs).  Therefore  we  can  apply  Theorem  3.2  to  establish  the  existence, 
uniqueness  and  continuous  dependence  of  a  global  solution  for  the  TCE  model  given 
by  (2)  -  (5). 

4  Estimation  of  parameters  and  convergence  of 
Galerkin  approximations 

The  TCE  PBPK-hybrid  model  is  dependent  on  many  physical,  biological  and  chem¬ 
ical  parameters,  some  of  which  can  be  explicitly  measured.  However,  most  of  the 
adipose  model  parameters  are  unknown,  and  it  would  be  impractical  to  measure 
them  experimentally.  Therefore,  we  must  use  inverse  problem  techniques  to  esti¬ 
mate  these  parameters. 

In  this  section  we  address  theoretical  issues  related  to  the  standard  parameter 
estimation  problem  associated  with  the  TCE  PBPK-hybrid  model,  and  we  establish 
the  theoretical  convergence  of  our  numerical  approximation  scheme.  The  results  we 
present  here  are  based  on  the  work  of  Banks  and  Kunisch  [4] ,  who  developed  results 
for  a  general  parameter  estimation  problem.  Banks  and  Musante  [6]  applied  these 
results  to  the  abstract  class  of  nonlinear  parabolic  systems  discussed  in  [5] .  Here  we 
extend  the  ideas  presented  in  [4]  and  [6]  to  the  TCE  model  and  the  abstract  class 
of  nonlinear  parabolic  systems  discussed  in  Section  3. 

Consider  the  abstract  problem 

.'/(/)  +  A(q)y(t)  +  g(q)(y(t))  =  f{t;  q)  (40) 

y(  0)  =  Vo,  (41) 

which  is  a  parameterized  formulation  of  the  system  (7),  (8).  Specifically,  the 
sesquilinear  form  a,  the  nonlinearity  g  and  the  forcing  function  /  are  now  func¬ 
tions  of  a  vector  q  of  parameters  which  must  be  estimated  using  data. 

The  parameter  vector  q  belongs  to  a  set  Q  of  admissible  parameters,  where  Q 
may  be  an  infinite  dimensional  space.  In  the  case  of  our  TCE  model,  the  physiology 
of  adipose  tissue  (see  [2,  22])  suggests  that  many  of  the  adipose  model  parameters 
may  be  dependent  on  time  and/or  space.  This  includes  parameters  such  as  the 
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volumetric  blood  flow  rate  Qf  to  adipose  tissue  and  the  permeability  coefficients 
H'BA,  Ubi  and  /ija-  Such  a  dependency  on  time  and/or  space  would  lead  to  infinite 
dimensional  parameters  and  an  infinite  dimensional  parameter  space  Q. 

Using  arguments  and  modified  conditions  as  outlined  in  Section  3,  we  can  readily 
establish  the  existence  of  a  unique  solution  to  the  system 

=  ( /(*;/?), $)v*,v  (42) 

2/(0)  =  yo  (43) 

for  ib  G  V.  Specifically,  the  modified  conditions  used  to  prove  the  existence  of  a 
solution  to  the  parameterized  system  (42),  (43)  include  the  extra  requirement  of 
uniformity  with  respect  to  the  parameters  q. 

These  solutions  y{-,q)  to  the  parameterized  problem  can  be  used  with  experi¬ 
mental  data  to  estimate  the  parameters  q  by  solving  the  least  squares  parameter 
estimation  problem  given  by 

Nt 

min  J(q,z)  =  V  | Oy(tit  •;<?)  -  zt \2  ,  (44) 

where  Zj.  i  1 ,  A/  are  observations  taken  at  time  t;,  and  | •  |  is  an  appropriately 
chosen  Euclidean  norm.  The  solutions  y(tj,-;q)  belong  to  an  infinite  dimensional 
state  space  T-L.  As  discussed  in  [4],  the  observation  operator  O  is  determined  by 
the  type  of  data  that  is  being  collected.  Examples  appropriate  for  the  TCE  model 
and  related  experiments  include  observations  of  concentrations  taken  in  time  at  a 
given  point,  or  observations  in  time  of  an  average  concentration  over  a  region.  This 
latter  type  of  observation  may  be  used  when  concentrations  are  measured  from 
homogenized  tissue  samples. 

As  in  [4,  6,  7],  we  address  the  parameter  estimation  problem  (44)  by  study¬ 
ing  the  convergence  properties  of  the  finite  dimensional  approximating  parameter 
estimation  problems 

Nt 

min  JA  (q)  =  ^  |  OyN  (U ,  • ;  q)  -  z,  \  ~  ,  (45) 

qeQM  “  1  1 

1  =  1 

where  the  sets  QM  are  a  sequence  of  finite  dimensional  sets  that  approximate  Q,  and 
j/A  (i)  G  %A  C  %  is  the  solution  to  the  following  finite  dimensional  approximation 
of  (42),  (43): 

(yN{t),ip)  +  (A{q}yN{t),ip)v*,v  +  (g{q){yN =  </(*;  q),4)v*,v  (46) 

yN(  0)  =  PNVo  (47) 

for  ijj  G  T-LN .  The  set  7iN  is  a  finite  dimensional  subspace  of  the  state  space  T-L,  and 
PA  is  the  orthogonal  projection  of  TiN  onto  Ti. 
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Assuming  certain  compactness  conditions  that  we  detail  below,  we  are  guaran¬ 
teed  that  for  each  M  and  N  there  exists  a  minimizer  qM’N  for  the  finite  dimensional 
estimation  problem  (45).  In  this  section  we  address  the  convergence  properties  of 
this  sequence  of  minimizers  as  M,  N  — >  oo,  and  the  circumstances  for  which  we  can 
guarantee  convergence  to  a  minimizer  for  the  infinite  dimensional  parameter  esti¬ 
mation  problem  (44).  The  results  we  prove  below  for  the  abstract  class  of  systems  in 
Section  3  are  adaptations  and  extensions  of  results  for  the  general  estimation  prob¬ 
lem  discussed  in  [4]  and  the  estimation  problem  for  the  class  of  nonlinear  parabolic 
systems  presented  in  [6] .  We  also  show  that  the  results  established  here  apply  specif¬ 
ically  to  the  TCE  PBPK-hybrid  model,  which  is  a  member  of  the  class  of  systems 
described  above. 

We  obtain  as  a  special  case  (with  qA  =  q  fixed  in  Theorem  4.3  below)  the 
convergence  of  solutions  for  the  finite  dimensional  approximations  (46),  (47)  to 
the  solution  for  the  infinite  dimensional  system  (42),  (43).  This  establishes  the 
convergence  of  the  numerical  scheme  that  we  implement  computationally  in  [22]  to 
generate  model  simulations  and  for  use  in  parameter  estimation  problems. 

4.1  Well-posedness  of  solutions  for  the  parameterized  system 

In  this  section  we  address  the  well-posedness  of  solutions  for  (40),  (41).  That  is, 
using  the  arguments  detailed  in  Section  3,  we  may  prove  for  this  abstract  system 
the  existence  of  a  unique  solution  which  depends  continuously  on  the  initial  data. 

As  in  Section  3.1,  we  define  the  state  spaces  V  and  'H  which  form  the  Gelfand 
triple 

where  V  is  continuously  and  densely  embedded  in  'H.  Moreover,  there  exists  k  >  0 
such  that 

\^\<mv  w 

for  all  ^  £  V,  where  |-|v  is  the  norm  on  V  and  |-|  denotes  the  norm  on  'H.  We 
continue  to  use  the  notation  (v)v*,v  to  represent  the  duality  product,  which  is  the 
extension  by  continuity  of  the  inner  product  (•,  •)  in  V.  (see  Section  3.1). 

Here  we  consider  the  system 

{y(t)$)v,v  +  p(q)(v(t),il>)  +  {g(q)(y(t)),il>)  =  {f(t;q)iip)v*,v  (49) 

1/(0)  =  y0  (50) 

for  ip  £  V,  where  the  sesquilinear  form  a(q)  :  V  x  V  — Z  R  is  associated  with  the 
operator  A(q)  :  V  — >  V*  so  that  cr(q)(u,v)  =  (A(q)u,  w)y,v  f°r  all  u,V  £  V. 
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We  make  the  following  assumptions  for  the  parameterized  system,  which  are 
the  same  assumptions  as  those  made  in  Section  3.1  with  the  extra  requirement  of 
uniformity  in  q  £  Q: 

(Al')  The  sesquilinear  form  a(q)  is  uniformly  bounded  in  V,  i.e.,  there  exists  C\  >  0 
independent  of  q  £  Q  such  that 

\v(q)(u,vj$;<  Ci  |u|v  |y|v  for  all  u, v  £  V.  (51) 

(A2')  The  sesquilinear  form  a (q)  is  uniformly  strictly  coercive  on  V.  That  is,  there 
exist  k\ ,  A  >  0  independent  of  q  £  Q  such  that 

a (q)(u,u)  >  k\  |w|y  —  A  |w|"  for  all  u  £  V.  (52) 

(A3')  The  forcing  function  f(-\q)  satisfies 

/(■;</)  €L2((0,T);V*)  (53) 

for  each  q  £  Q. 

(AT)  The  nonlinear  function  g(q)  :  TT  — >  Ti  satisfies  the  following  uniform  local 
Lipschitz  condition:  let  Br( 0)  =  {u  £  T-L  :  |u|  <  r}  denote  the  ball  of  radius 
r  centered  around  the  origin  in  T-L.  Then  given  r  >  0,  there  exists  Lb,.  >  0 
independent  of  q  £  Q  such  that 

| g(q)(u)  -  g(q)(v)\  <  LBr  \ u  -  u|  for  all  u,v  £  Br(0).  (54) 

(A5')  There  exist  positive  constants  CA  and  C3  independent  of  q  £  Q  such  that 

\g(q)(u)\  <  Co  |u|  A  C3  for  all  u  £  TL.  (55) 

Under  these  assumptions  we  can  prove  the  existence  of  a  unique  solution  to  the 
parameterized  system  (49),  (50).  The  proof  of  this  theorem  is  similar  to  the  proof 
for  Theorem  3.2. 

Theorem  4.1  Under  Assumptions  (Al1)  -  (A51)  and  for  any  yo  £  Ti,  the  sys¬ 
tem  (49)  -  (50)  has  a  unique  global  solution  y  £  Lo((0,  T);  V)  fl  £7( [0,  T}\ TL)  which 
depends  continuously  on  the  initial  data ,  with  y  £  I/o((0,  T);  V*). 

4.2  The  general  parameter  estimation  problem 

In  this  section  we  present  theorems  related  to  convergence  for  the  general  parameter 
estimation  problem 

Nt 

min  J(q,z)  =  V  | Oy(tis  -,q)  -  Zi\2  , 
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given  the  observations  Zj,i  =  1, . . . ,  Nt  taken  at  time  t,  and  the  observation  operator 
O.  Moreover,  •;  q)  are  the  solutions  of  the  parameterized  system  (49),  (50).  As 
in  Section  3.1,  we  use  Galerkin  type  approximations  applied  to  the  parameterized 
system  (40),  (41). 

Proceeding  as  in  [6,  7],  we  make  the  following  assumptions  about  the  state  and 
parameter  spaces  TL,  Q,  TLN  and  QM: 

(Bl)  The  sets  Q  and  QM  each  lie  in  a  metric  space  Q  with  metric  d.  We  assume 
that  Q  and  QM  are  compact  in  this  metric,  and  that  there  is  a  mapping 
iM  :  Q  — >  QM  such  that  QM  =  iM(Q).  Moreover,  for  each  q  G  Q  we  have 
iM{q )  — >  q  in  Q,  where  the  convergence  is  uniform  in  q  G  Q. 

(B2)  The  finite  dimensional  subspaces  Tt A  satisfy  TL N  C  V  for  all  N. 

(B3)  For  each  if  G  V,  | ip  —  PN tjj |  — y  0  as  N  — >  oo,  where  PA  is  the  orthogonal 
projection  operator  of  %  onto  H.N . 

Moreover,  we  assume  that  A(q),  g(q)  and  f(-;q)  depend  continuously  on  the  pa¬ 
rameters  q  G  Q  in  the  sense  that 

(Cl)  |<7 (q)(u,v)  —  a(q)(u,.v)\  <  d\{q,q)  |u|v  |n|>p  for  all  u,v  G  V,  where  d\(q,q)  — >  0 
as  d(q,  q)  0. 

(C2)  \g(q)(u)  -  g(q)(u)\  <  do(q,q )  |w|  for  all  u  G  P,  where  d2(q,q)  ->  0  as  d(q,q )  -> 
0. 

(C3)  The  mapping  q  — y  f{-;q )  is  continuous  from  Q  to  I/o((0,  T);  V*). 

Note  that  the  compactness  conditions  of  Assumption  (Bl)  guarantee  the  exis¬ 
tence  of  a  minimizer  qM’N  to  the  finite  dimensional  estimation  problem  (45).  The 
following  theorem  from  [7]  (Theorem  5.1,  p.  124)  establishes  conditions  for  which 
the  sequence  of  minimizers  qM’N  converges  to  a  minimizer  for  the  infinite  dimen¬ 
sional  estimation  problem  (44). 

Theorem  4.2  To  obtain  convergence  of  at  least  a  subsequence  of  {qM,N}  to  a  solu¬ 
tion  q  of  minimizing  (44)  subject  to  (4%),  (43),  it  suffices,  under  Assumption  (Bl), 
to  argue  that  for  arbitrary  sequences  {qM,N}  in  Q  with  qM’N  — >  q  in  Q,  we  have 

OyN(t;qM’N)^Oy(t;q).  (56) 

The  condition  (56)  implies  that  the  original  sequence  of  Galerkin  approxima¬ 
tions  in  Section  3.1  converges,  which  is  stronger  than  the  subsequential  convergence 
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proven  there.  This  stronger  result  is  important  for  establishing  the  convergence  of 
the  numerical  scheme  that  we  develop  and  implement  in  [22], 


Now  we  verify  the  conditions  of  Theorem  4.2  for  the  general  estimation  problem 
described  in  this  section. 

Theorem  4.3  Assume  that  Assumptions  (Bl)  -  (B3)  and  (Cl)  -  (C3)  are  satis¬ 
fied,  in  addition  to  the  conditions  of  Theorem  f.l.  Let  be  any  sequence  in  QA 
satisfying  <jA  — >■  q  G  Q.  Then 

yN (t;  qN)  — >  y(t;  q)  in  TL  uniformly  on  [0,T], 

yA  (t;  <jA  )  — >  y(t ;  q)  in  V  for  almost  all  t  >  0, 
where  yN  satisfies 

(yN(t),C)  +  &(qN)(yN  (t),if)  +  (cj(qN)(yN (t)),if)  =  </(/;  r/ Vh r)v,v  (57) 
yN(  0)  =  PNyo  (58) 

for  all  if  E  TLN ,  arid  y  satisfies  (49)  -  (50)  for  all  if  €  V. 

Proof:  It  follows  from  Theorem  4.1  that  y(t)  &  TL  for  every  t  G  [0,  T]  and  y(t)  G  V 
for  almost  every  t  G  [0,T].  Moreover, 

</A7)  -  y(t;  g)|  <  \yN(t;  qN )  -  g)|  +  \PNy(t;  q)  -  y(t ;  g)|  . 

By  Assumption  (B3),  we  have  —  PN  if  |  — >  0  for  each  if  G  V  as  Ar  — >■  oo. 

Note  that  y  G  C([0,T];P)  implies  that  {y{t)  :  t  G  [0,T]}  is  compact  in  TL,  which 
along  with  (48)  and  the  dense  embedding  V  ^  TL  further  implies  that 

i,x!j(T-<i)  ~y(t;q)\  >  0 

uniformly  in  t  G  [0,T],  Therefore  it  suffices  to  show  that 

\yN(t',qN)  —  PNy(t',q)\  — >  0 
uniformly  on  [0,T]  as  N  — >  oo. 

As  in  [6],  we  define 

yN  =  yN (t;  qN) 

y  =  y(t-,  q) 

AN  =  y»(t:qN)-PNy(t:q)=  yN  -PNy. 

It  follows  that 

A"  =  i,K  -  ftP"y, 
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and  hence  for  all  ip  G  Ji A  we  have 

(A N,iP)v*,v  =  {yN  -  ^P*y,tp)v*,v 

=  (yN  -y  +  y-  ^j.pNy^  V>)v*,v 
=  -  {y,tp)v,v  +  (y  -  ^P*y,tp)v*,v- 

Moreover,  we  use  (49)  and  (57)  to  obtain 

(A N,ip)v*,v  =  {y  -  ^PNy,iP)v*,v  +  if(t',qN)  -  f(t\q),tp)v*,v 
+  o{qN){y  -  PNyi*P)  +  v(q)(y,ip)  -  o{qN){y,ip) 

+  {g(q)(y)-g(qN)(yNU)  -  °(qN)(yN  -  PNy ,4) 

for  all  ip  G  JiN ,  where  we  have  added  and  subtracted  cf(qN)(y  —  PNy,ip).  It  then 
follows  from  the  definition  of  AA  that 

(A N ,ip)V\v  +  (?(qN)(AN  ,ip)  =  {y  -  ^pNy,^)v*,v  +  {f{t\qN )  -  /(*;«), ^)v*,v 
+  <j(q)(y,  ip)  -  a(qN)(y,ip)  4-  a {qN)(y  -  PNy,ip) 

+  (g(q)(y)-g(qN)(yN),i’)  (59) 

for  all  ip  G  PA .  Now  we  set  ip  =  AA  G  PN  in  (59)  to  arrive  at  the  equation 

\jt  [|4AT]  +  <r(qK)(&K,&N)  =  (y  -  lP"y,  A>.,v 

+  ( f(t;qN )  -  /(i;g),  AJV)v.,v  +cr(g)(y,  Aw)  -  a(qN)(y,  AN) 

+  ^(gAr)(y  -  J^J/,  AAr)  +  <5(g)(j/)  -  g(qN)(yN),  AN)  (60) 

for  almost  all  t  G  [0,T].  We  integrate  the  left  side  of  (60)  and  apply  Assumption 
(A2')  and  the  initial  condition 

AAr(0)  =  /r(0)  -  PNy( 0)  =  yN( 0)  -  PNy0  =  0 

to  obtain  the  inequality 

J*  ~\&N(s)f  ds  +  J*  a(qN)(AN(s),AN(s))ds 

>  ^  |Aa  (t)\~  +  k\  I  | AiV(s) | ~  ds  —  A  f  |AA  (s)|“  c(ill) 

2  Jo  Jo 

for  all  t  G  [0,T]. 

Next  we  integrate  the  right  side  of  (60)  from  0  to  t.  Note  that  the  first  term 
is  equal  to  zero,  which  follows  from  the  fact  that  PA  is  the  orthogonal  projection 
from  7i  onto  .  Indeed,  since  y  —  PA  y  G  Ji  we  have 

(y  -  J \pNy > ^)v*,v  =  ^.{y  ~  P A !)■. i’)  =  0  (62) 
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for  all  tp  £  since  PA  is  the  orthogonal  projection  of  H  onto  7iN ,  which  implies 
that  (I  —  PN)y  is  orthogonal  to  all  elements  in  7 iN  (including  A N(s)  for  all  s  £ 

[0  ,T}). 


Integrating  the  second  term  on  the  right  side  of  (60),  we  use  the  inequality 
ab  <  (1/2 )ar  +  (1/2 )b'2  to  obtain 

J  (f(s-,qN)  -  f{s-,q),AN{s))v,,vds  <  J  \f(s\qN)  -  f(s]  q)\~v,  ds 

+  ^  I AiV(s) | v  (63) 

where  j3  >  0.  We  proceed  in  a  similar  manner  with  the  third  and  fourth  terms  on 
the  right  side  of  (60),  applying  Assumption  (Cl): 

rt  pi  N\  ft 

^  [v(q)(y(s),  AA  (4)  -  o{qN){y{s),  AN (s))]  ds  <  - j ^  [y(s)\2v  ds 

+  |  |AAr(s)|y  ds-  (64) 

The  integral  of  the  fifth  term  on  the  right  side  of  (60)  can  be  bounded  using  As¬ 
sumption  (Al'): 

<J(qN)(y(s)  -  P* y(s),  AN  (s))ds  <  ^  ^  | y(s)  -  PNy(s)\~v  ds 

+  \  |AA  (s)|y  ds.  (65) 

Finally,  we  integrate  the  sixth  term  from  0  to  t,  obtaining 

[\g(q)(y(s))-  g(qN)(yN(s)),AN(s))ds 
Jo 

<  f  \(g(q)(y(s))  -  g(qN)(y(s)),AN(s))\ds 
Jo 

+  f  \(g(qN)(y(s))  - g(qN)(yN(s)),AN(s))\ds.  (66) 

J  o 

The  first  term  on  the  right  side  of  (66)  can  be  estimated  using  Assumption  (C2) 
and  (48): 


\{g(q)(y(s))  -g(qN)(y(s)),AN(s))  ds  <  /  d2(q,qN)k2  \y(s)\ 


We  obtain  a  bound  for  the  second  term  on  the  right  side  of  (66)  using  Assumption 
(A4')  applied  to  y(s)  and  PN,y(s).  In  the  proof  of  Theorem  3.2,  we  showed  that  the 
solution  y  £  T2((0,  T);  V)  fl  C([0,T];'H)  to  the  system  (15),  (16)  is  bounded  by  the 
constant  C  as  in  (37).  It  is  easily  seen  that  for  each  s  £  [0,T],  the  same  bound 


\y(s)\2<c 
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applies  for  the  solution  y  £  L2((0,T);V)  fl  C([0, T\,T-L)  to  the  system  (49),  (50) 
which  is  guaranteed  by  Theorem  4.1.  Moreover,  for  each  s  £  [0,T]  we  have 


\PNy(s)\  <  \\PN\\  |j/(s)|  =  [y(s)\ , 

which  implies  that  both  y(s)  and  PNy(s )  are  in  the  ball  B  in  7~L,  Using  the  same 
arguments  and  the  bound  (36),  we  can  show  that  yN{$)  is  in  B  ^j-t.  Hence  we  may 
apply  Assumption  (A4')  and  (48)  to  obtain 

f  MqN){y{s))-g{qN){yN(s)),\N{s))\ds 


< 


+ 


< 


\(9(qN)(y(s))-g(qN)(PNy(s)),AN(s))\ds 
\(g(qN)(PNy(s))  -  g(qN)(yN(s)),  AN(s))\ ds 


k2L% 


Vs 


2/1 


| 2/00  -  P A y(s)\~v  ds  +  I  I AA" (a)|“  ds 


+  Lb 


Vs 


AA  ( s )  “  ds. 


(68) 


Now  we  combine  the  terms  in  (61)  -  (68)  and  use  (60)  to  arrive  at  the  estimate 
\  |AA(i)f  +  [k,  -  (3  +  2/A;i)  J* \AN(s)\2vds 


< 


+ 


+ 


j-j  Jo  |/(a;  gA)  -  f(s;  g)|y.  ds  +  d  "!2'l  ^  ^  \y{s%ds 
C'l  f4  |  N  ,2  ,  k2d2(q,  qN)  [4  o 

2d  Jo  “  p  y ^  I V  ds  + - ^ - Jo  Ms) I v  ds 


k2L% 


Vs 


2/1 


\y(s)  -  PNy(s%  ds  +  ( Lb ^  +x)  Jo  |A-(,)|-  ds 


=  SN(t)  +  (LBv_  +  X)  j  |AA  (s)|2  ds 
<  Sn(T)  +  (LBv_  +  \)  J*\AN(s)\2ds 


(69) 


for  t  G  [0,T],  Note  that  as  N  — »  oo  we  have  SN(t)  — *  0  for  each  t  £  [0,T].  This  is  a 
consequence  of  Assumption  (C3),  q A  — >■  q,  y  G  L2((0,T);  V),  and  Assumption  (B3). 


Now  we  choose  /i  >  0  in  (69)  so  that  k\  —  (3  +  2Ar)f  >  0.  We  then  apply 
Gronwall’s  inequality  to  obtain 

| AiV(t)|2  <  2e2(LBv^+A)T(5Ar(T)  ->■  0,  (70) 

which  implies  that  AN  — >  0  in  C([0,T];'H)  as  N  — >  oo.  Moreover,  combining  (69) 
with  (70)  we  have 

(ki  -(3  +  2fc2)0  J*\AN(s)\lds  <  6N(T)  +  (LBv_  +  \)  J*\AN(s)\2ds 
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<  Sn(T)  ^1  +  2T(Lb +  A)e2(iBVtf +A)Tj 
->■  0 

and  hence  AA  — >  0  in  Ln((0,T);  V).  This  completes  the  proof  of  the  theorem. 

4.3  Parameter  estimation  for  the  TCE  model 

The  TCE  PBPK-hybrid 

model  outlined  in  Section  2  includes  several  physical  parameters,  each  of  which  is 
considered  a  positive  constant  that  is  bounded  above  by  some  maximum  admissible 
value.  This  implies  that  the  parameters  for  this  model  lie  in  a  compact  subset  of 
Euclidean  space.  Therefore  we  set  QM  =  Q  for  all  M. 

Using  arguments  detailed  in  [22],  it  can  be  shown  that  Assumptions  (Al')  - 
(A5')  are  satisfied  for  the  TCE  PBPK-hybrid  model.  We  state  without  proof  that 
Assumptions  (Cl)  -  (C3)  also  are  satisfied  for  the  TCE  model.  Moreover,  it  is 
seen  that  Assumptions  (Bl)  -  (B3)  hold  for  our  choice  of  state  spaces  V,  77  and  7/A 
defined  in  Section  2.1.  It  follows  that  the  theoretical  results  presented  in  Section  4.2 
apply  to  the  TCE  PBPK-hybrid  model. 

Experimental  observations  collected  in  the  experiments  of  Evans  et  al.  [15]  in¬ 
clude  measurements  of  TCE  concentrations  in  the  venous  blood,  as  well  as  unbound 
concentrations  in  the  fat,  liver  and  brain  tissue.  The  tissue  concentrations  are  col¬ 
lected  from  homogenized  tissue  samples,  leading  to  a  measurement  of  an  effective 
average  tissue  concentration.  For  the  liver  and  brain  tissues  our  model  assumes 
well-mixed  compartments,  which  implies  that  the  average  tissue  concentration  is 
equal  to  the  concentration  at  any  point  in  that  tissue.  Therefore  the  observation 
operators  for  the  venous  blood  and  the  brain  and  liver  tissues  are  given  by 

Oy{ti,-,q)  =Cj(tnq),  (71) 

where  Cj  =  Cv,  Cbr  and  Ci  respectively. 

In  the  adipose  tissue  our  model  assumes  spatially-varying  concentrations,  so 
that  an  effective  average  concentration  is  calculated  by  taking  an  average  of  con¬ 
centrations  over  the  adipocyte  region.  This  implies  an  observation  operator  of  the 
form  ^ 

Oy{U,-,q)  =  j  CA{U,O,<j>)d0  d<f>.  (72) 

Each  of  these  observation  operators  is  utilized  in  [22],  where  we  present  detailed 
computational  parameter  estimation  results  for  the  TCE  model. 
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5  Concluding  remarks 


In  this  paper  we  have  addressed  theoretical  issues  related  to  the  TCE  PBPK-hybrid 
model  and  a  general  class  of  toxicokinetic  models.  Specifically,  we  established  ex¬ 
istence,  uniqueness  and  continuous  dependence  results  for  an  abstract  class  of  non¬ 
linear  parabolic  systems  which  includes  the  TCE  PBPK-hybrid  model  as  a  special 
case.  Moreover,  we  addressed  theoretical  issues  for  the  associated  general  parameter 
estimation  problem  and  established  the  convergence  of  the  Galerkin  finite  element 
approximations . 
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